Results from the First Science Run of the ZEPLIN-III Dark Matter Search 
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The ZEPLIN-III experiment in the Palmer Underground Laboratory at Boulby uses a 12 kg two- 
phase xenon time projection chamber to search for the weakly interacting massive particles (WIMPs) 
that may account for the dark matter of our Galaxy. The detector measures both scintillation and 
ionisation produced by radiation interacting in the liquid to differentiate between the nuclear recoils 
expected from WIMPs and the electron recoil background signals down to ~10 keV nuclear recoil 
energy. An analysis of 847kg-days of data acquired between February 27 th 2008 and May 20 th 2008 
has excluded a WIMP-nucleon elastic scattering spin- independent cross-section above 8.1 x 10 -8 pb 
at 60 GeVc -2 with a 90% confidence limit. It has also demonstrated that the two-phase xenon 
technique is capable of better discrimination between electron and nuclear recoils at low-energy 
than previously achieved by other xenon-based experiments. 
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A. Motivation 

Searches for weakly interacting massive particles 
(WIMPs) are motivated by the coming together of unifi- 
cation schemes, such as supersymmetry, which predict 
new particle species, and extensive observational evi- 
dence which demonstrates the need for additional non- 
baryonic gravitational mass within the Universe. That 
the WIMPs of supersymmetry naturally fulfill this need 
is remarkably persuasive. Indeed, WIMPs occur in other 
frameworks too. As a generic class of particle they are 
assumed to only interact non-gravitationally with bary- 
onic matter via the weak interaction. Whilst this offers 
a mechanism for energy transfer and hence detection, it 
also implies rather low event rates and energy deposits: 
<0.1 events/day/kg and <50 keV respectively. This dic- 
tates the use of sensitive underground experiments capa- 
ble of specifically identifying energy deposits due to elas- 
tic scattering of incoming particles from target nuclei. 
ZEPLIN-III is the latest in a progressive series of instru- 
ments designed to push steadily the sensitivity limits by 
exploring alternative approaches using xenon-based tar- 
gets [1, 2]. 
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ZEPLIN-III is a two-phase (liquid/gas) xenon time- 
projection chamber specifically designed to search for 
dark matter WIMPs. Its design and performance details 
have already been presented elsewhere [3, 4] and only a 
brief reminder is given here. The experiment is operat- 
ing 1100 m underground. The active volume is a disc of 
35 mm thickness and ~190 mm diameter which contains 
~12 kg of liquid xenon above an array of 31 2-inch di- 
ameter photomultipliers (PMTs). The PMTs employed 
during this first science run were ETL D730/9829Q [5], 
and they were used to record both the rapid scintillation 
signal, SI, and a delayed second signal, S2, produced by 
proportional electroluminescence in the gas phase above 
the liquid [6] . The PMT array was immersed in the liquid 
viewing upwards. The electric field in the target volume 
was defined by a cathode wire grid 36 mm below the liq- 
uid surface and an anode plate 4 mm above the surface 
in the gas phase. These two electrodes alone produce the 
drift field in the liquid (3.9 kV/cm), the field for extrac- 
tion of the charge from the surface, and the electrolumi- 
nescence field in the gas (7.8 kV/cm). A fiducial volume 
for WIMP searches was defined by using a time window 
for delays between SI and S2, which selected a depth 
slice within the liquid, and by 2-D position reconstruc- 
tion from the PMT signals to select a radial boundary 
at 150 mm. The time window was set between 500 ns 
and 13,000 ns which selected depths between 1.29 mm 
and 33.43 mm. These together defined a fiducial volume 
containing 6.5 kg of xenon. 
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The PMT signals were digitised at 2 ns sampling over a 
time segment of ±18 /is either side of the trigger point. 
Each PMT signal was fed into two 8-bit digitisers (AC- 
QIPJS DC265) with a x 10 gain difference between them 
provided by fast amplifiers (Phillips Scientific 770), to 
obtain both high and low sensitivity read-out covering 
a wide dynamic range. The PMT array was operated 
from a common HV supply with attenuators (Phillips 
Scientific 804) used to normalise their individual gains. 
The trigger was created from the shaped sum signal of 
all the PMTs. For nuclear recoil interactions the trig- 
ger was always caused by an S2 signal for energies up 
to Sl=40 keVee, where keVee is an energy unit refer- 
enced to the equivalent SI signal produced by 122 keV 
7-rays from 57 Co. The trigger threshold was ~11 ioni- 
sation electrons and this corresponded to ^0.2 keV for 
electron recoils (for nuclear recoils see Section HID 2). 
This S2 threshold was set to avoid excessive triggers from 
single electron emission events and from electron and nu- 
clear recoils whose primaries would otherwise have been 
undetectable as they fall below the SI detection thresh- 
old. 

The xenon target was contained within a vessel itself 
located within a vacuum jacket both made from low- 
background oxygen-free copper. Cooling was provided 
by a 40 litre liquid nitrogen reservoir, also made from 
copper, inside the vacuum jacket. Thermal stability to 
<0.5 °C was achieved over the entire run by controlling 
the flow of cold nitrogen boil-off gas through the base- 
flange of the xenon vessel. Pressure stability to 2% was 
maintained. The ZEPLIN-III detector was completely 
surrounded by a shield of 30 cm thick polypropylene and 
20 cm thick lead, giving 10 5 attenuation factors for both 
7-rays and neutrons from the cavern walls. Dedicated ac- 
cess through the shield was provided for the radioactive 
calibration source delivery, instrument levelling screws 
and pipe-work to the external gas purification system. 



C. Science Data 

WIMP-search data were collected over 83 days of 
continuous operation in the Boulby Laboratory starting 
on 27 th February 2008. An 84% live time was achieved 
during the science run and some 847kg-days of raw 
data were collected from the 12 kg target volume. 57 Co 
calibration measurements were made every day. Nuclear 
recoil calibrations were made with an AmBe neutron 
source at the beginning and end of the 83 day period 
(5 hrs each). A typical event, from a neutron elastic 
scattering interaction in the liquid with Sl=5 keVee, 
is shown in Figure 1 as recorded through the high- 
sensitivity sum channel. A short Compton calibration 
was performed using a 137 Cs source at the beginning of 
the run with a much longer run at the end (122 hrs). 
Ten percent of the science data (every 10 th file) were 
used to develop initial data analysis and selection cuts, 
to establish the level of the electron-recoil background, 
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FIG. 1: Segment of the high-sensitivity summed waveform for 
a neutron elastic scattering event with Sl= 5 keVee, showing 
a small primary pulse (SI) preceding a large secondary pulse 
(S2). Some PMT after- pulsing and, possibly, single electron 
emission can be seen following S2. Note that only excursions 
>3 rms on individual channels are added into the summed 
waveform. See later text for more detailed discussion of some 
of these points. 



and to define the boundaries for the WIMP-search box 
and its acceptance. At first, the remaining 90% of the 
science data were retained unopened to carry out a 
'blind' analysis, but these data were eventually used for 
perfecting some data-selection cuts as detailed below, 
making the final analysis non-blind. 

Pulse-finding algorithms were used to identify 
signals in the 62 waveforms (independently for each 
PMT and for high and low sensitivity channels). These 
were then categorised as SI or S2 candidates based on 
a pulse width parameter (charge mean arrival time, r): 
scintillation pulses are much shorter (r<40 ns) than 
electroluminescence pulses, with durations correspond- 
ing to the drift time across the gas gap (r~550 ns). 
Viable SI and S2 candidates were then subject to 
software thresholds (>3 channels recording signals above 
1/3 photoelectron (p.e.) for SI and a minimum area of 
~5 ionisation electrons for S2). Only events with one 
SI and one S2 were considered for further analysis. Of 
particular note here, x 2 goodness of fit indicators within 
the position reconstruction of both SI and S2 were 
used to remove multiple-scatter events, and this was 
particularly effective for those with one vertex in a 'dead' 
region of the xenon, which would otherwise have been 
a troublesome background. Such 'dead' regions include 
the reverse-field volume between the cathode wire and 
the PMT grid wire [4] and the thin (0.5 mm) layer of 
xenon surrounding the PMT bodies. Double-Compton 
interactions with at least one vertex in these regions, 
referred to as 'multiple-scintillation single-ionisation' 
(MSSI) events, fulfil the previous selection criteria 
since there is no S2 pulse from the dead region and 
the coincident scintillation pulses are added together 
in a single SI. Unfortunately, perfecting this selection 
eventually required use of the full data-set as will be 
described in more detail below. 
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II. CALIBRATION 

A. Scintillation Response and Position 
Reconstruction 

An external 57 Co source was inserted through the 
shield and located above the instrument every day. The 
dominant 122 keV 7-rays have a photoelectric absorp- 
tion length of 3.3 mm in liquid Xe, and hence provided 
good standard calibration candles from interactions close 
to the liquid surface. A typical 57 Co spectrum is shown 
in Figure 2. The SI signal channel exhibited a light de- 
tection efficiency at our operating field (3.9 kV/cm) of 
L y =1.8 p.e./keVee, decreasing from 5.0 p.e./keVee on 
application of the electric field. The 122 keV interac- 
tions were used for a number of purposes to calibrate 
the instrument. Using S2 pulses, an iterative procedure, 
whereby a common cylindrical response profile was fit- 
ted to each channel, was used to normalise the mea- 
sured response from each PMT (i.e. 'flat-field' the ar- 
ray). Position reconstruction in the horizontal plane was 
then achieved by using the converged response profiles in 
a simultaneous least-squares minimisation to all chan- 
nels [7]. This method complements the Monte Carlo 
template matching procedure also being used but is less 
dependent on accurate iterative simulations [8]. Finally, 
the integrated areas of the SI and S2 responses gave light 
collection correction factors as a function of radial posi- 
tion. Using this procedure a full-volume energy resolu- 
tion of (7=5.4% at 122 keV was obtained with an energy 
reconstruction using a combination of the SI and S2 re- 
sponses to reflect the fact that, for electron recoils, these 
two channels are anti-correlated at a microscopic level. 
The individual SI and S2 resolutions at 122 keV are 
16.3% and 8.8%, respectively. Also shown in Figure 2 
is the comparison of the response to simulation. Not 
only are the two main 57 Co lines well fitted but there is 
also a good match to the predicted Compton feature at 
~35 keV. The excess above 150 keV is mainly due to the 
unsubtracted background. The left-hand panel in Fig- 
ure 3 shows the distribution in the x-y plane of events 
seen from the source. As expected most events are lo- 
cated towards the centre (the offset is due to an offset 
source position) with a radial fall-off as expected from 
the increasing thickness of copper along the line of sight. 



B. Stability, Electron Lifetime and Detector Tilt 

The 57 Co daily calibrations were used to assess the 
evolution of other operational parameters over the entire 
run: i) the average light and ionisation yields, as mea- 
sured by fits to the 57 Co SI and S2 pulse area spectra; 
ii) the mean electron lifetime in the liquid, obtained from 
the exponential depth dependence of the ratio of the ar- 
eas of the S2 and SI signals (hereafter simply referred 
to as S2/S1); iii) the evolution of the long-term detector 




energy (0.35*S1+0.65*S2), keVee 

FIG. 2: Response to an external 57 Co 7-ray source in the com- 
bined energy channel, exploiting SI and S2 anti-correlation. 
One day's experimental data are shown in blue with statis- 
tical error bars. The simulation result is indicated in red: 
the solid histogram shows the bare energy deposits and the 
shaded one shows the result of Gaussian-smearing with the 
energy resolution indicated in the figure. 
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FIG. 3: Distribution in the horizontal plane of events from 
the 57 Co source on the left and from the AmBe source on the 
right. The source positions are different for each image and 
neither is centred. In both cases the volume distribution is 
as expected from Monte Carlo simulations, given the location 
of each source. Interaction vertices can be seen out to the 
edge of the fiducial volume at a radius of 150 mm (red circle). 
The outer circle shows the edge of the liquid xenon target. 
Each PMT is marked by two smaller circles (PMT centres 
and envelopes). 



tilt due to local geological factors, as given by the polar 
dependence of the S2-width distribution, which probes 
the thickness of the gas layer. The detector tilted by 
less than 1 mrad over the run. Over the fiducial vol- 
ume this corresponds to a systematic change in the gas 
gap of < 3.5%, which in turn translates proportionally 
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into a variation in the S2 signal. This was not deemed 
sufficient to warrant a full correction [31]. The scintilla- 
tion mean light yield remained stable to a few percent, 
as did the ionisation yield, after correcting for the elec- 
tron lifetime in the liquid. Remarkably, the lifetime did 
show an evolution during the run in the form of an im- 
provement: from an initial value of 20 fj,s, achieved by 
initial gas-phase purification through external getters, a 
value of 35 /is had been reached by the end of the run 
(the full drift length of the chamber is only 14 fj,s) . There 
was no active recirculation used and this improvement is 
attributed to the clean, xenon-friendly materials used in 
detector construction and to the uninterrupted applica- 
tion of the electric fields during the entire run. As the 
area ratio S2/S1 is the main discriminant between nu- 
clear and electron recoils, a depth-dependent correction 
must be applied to the S2 area to compensate for electron 
trapping by impurities. The electrons from the deepest 
events within the fiducial volume drifted for 13 fj,s and the 
correction factor for these varied from 1.92 at the start of 
the run to 1.45 at the end. The daily 57 Co calibrations 
allowed this to be monitored throughout the science run 
and events were corrected individually using an historical 
trend profile. 



C. Linearity 
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FIG. 4: Expected mean number of SI photo-electrons as a 
function of the mean pulse area observed in the central chan- 
nel in the array. The expected signal is the mean of the Pois- 
son distribution obtained by counting the frequency of 'zeros', 
i.e. the absence of any response. 



The linearity of the response of each channel in the ar- 
ray was investigated using low-energy Compton-scattered 
events from the 137 Cs source, in order to rule out hard- 
ware and software distortion for processing of small sig- 
nals. The position of the vertex for each interaction was 
found and the waveforms from PMTs located a certain 
distance away from the vertex were selected based on the 
expected number of SI photons, given the cylindrical re- 
sponse profile determined from the 57 Co data as pointed 
out in II A. Provided that the expected number is in- 
deed small, the mean of the Poisson distribution for the 
number of detected photons can be quite accurately de- 
termined by counting the fraction of waveforms which 
do not contain any identified pulses, i.e. the frequency 
characterising the absence of any signal. This assertion 
is made against a sample of pure noise in the same wave- 
form. Repeating this procedure for all channels and a 
range of expected signal allowed comparison of the mean 
SI pulse area recorded in each trial against the expected 
Poisson mean, as shown in Figure 4 for the central PMT. 
In addition, this provides a very robust method to ob- 
tain the mean size of one photoelectron [9]. This has 
been calculated for every PMT within the array: the re- 
lationship is found to be linear to within the statistical 
accuracy of the measurement over a factor of 10 in mean 
pulse area, which covers the range of interest for WIMP 
nuclear recoil signals. The slope of the line in Figure 4 
provides a measure of the mean single photoelectron re- 
sponse (SER) for that PMT. The mean SER of all the 
PMTs in the array has been found in this way to be in 



the range 47 ± 12 pVs. The spread in these values forms 
part of the 'flat-field' correction discussed earlier; other 
dominant factors are the PMT quantum efficiency and 
imperfect hardware equalisation. 

D. Nuclear Recoil Response 

The nuclear recoil response in the energy range of in- 
terest to WIMP signals has been calibrated with neutrons 
from an AmBe (a,n) source. The source was placed in- 
side the polypropylene shielding above the detector but 
displaced to one side to reduce the interaction rate. The 
right-hand panel in Figure 3 shows the reconstructed 
event positions from the second calibration performed 
just after the science run had been completed. The dis- 
tribution is slightly non-unform in the x-y plane as ex- 
pected. 

Figure 5 shows a 'scatter-plot' of \og w (S2/Sl) as a 
function of energy in keVee from the AmBe calibration. 
The red line shows a smooth fit to the median of the 
elastic scatter distribution with ±lcr boundaries as blue 
lines. To obtain these curves the data were histogramed 
into 1 keVee bins and fitted by log-normal distributions. 
Examples of the quality of the fits are shown in Figure 6. 
The other well defined population in Figure 5, between 
40-70 keVee, is due to inelastic scattering of neutrons 
from 129 Xe nuclei and the more diffuse horizontal pop- 
ulation is caused by associated 7-ray interactions. The 
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FIG. 5: Calibration of the nuclear recoil response with an 
AmBe neutron source, plotted as the discrimination param- 
eter (log 10 (52/51)) as a function of 'electron-equivalent en- 
ergy' (i.e. using the SI channel calibrated by 57 Co). The 
lines show the trends of the mean and standard deviation of 
energy-binned log-normal fits to the recoil population. The 
distinct population above ~40 keVee is due to inelastic neu- 
tron scattering off 129 Xe nuclei. 

elastic nuclear recoil median turns out to be very closely 
approximated by a power law, which is shown most effec- 
tively by replotting the figure in log-log form (Figure 7). 
Not only is the power-law behaviour very apparent but it 
can also be seen that there is less obvious flaring at lower 
energies than seen in other xenon experiments whose data 
were taken at much lower electric fields [2, 10]. Also 
shown are lines illustrating the approximate thresholds 
for SI and S2. 



E. Electron Recoil Response 

The electron recoil response at low energies was 
established using a long duration calibration with a 
137 Cs radioactive source. Compton scattering of the 
662 keV 7-rays produced a significant number of events 
down to ~ 2 keVee but with only a small number 
extending far enough down in the S2/S1 parameter to 
reach the nuclear recoil median (Figure 8). The general 
behaviour of the electron recoil band is reminiscent of 
the XENON10 results [10-12], but with a slightly more 
pronounced upturn at low energy, a larger separation 
between electron and nuclear recoil bands and narrower 
distributions. The low-energy electron-recoil populations 
in the 137 Cs and the WIMP-search data-sets were fitted 
in 1 keVee bins by a skew-Gaussian function. The 



fits were performed using a maximum likelihood (ML) 
method with a Poisson distribution as estimator for the 
observed data. Three of the fits are shown in Figure 6. 
The distribution parameters are consistent bin-by-bin 
for the 137 Cs and WIMP data-sets, as confirmed in 
Figure 8. However, there are two distinct differences in 
the general behaviour. Firstly, the mean of the 137 Cs 
data is systematically lower than that of the WIMP 
data. It has been shown that this reduction is due to the 
high count rate used in collecting the Cs data causing 
the gain of the PMTs to be slightly suppressed due 
to saturation effects at low-temperature [9]. However, 
it was not feasible to lower the rate and still acquire 
sufficient data in a reasonable time and uncontaminated 
by other background. Secondly, the behaviour of the 
137 Cs data-set in the low S2/S1 tails is not closely 
representative of the science data, with the former 
exhibiting significantly more outliers. These events are 
attributed to MSSI double-Compton events as had been 
anticipated in [3]. This is not evident in Figure 8 as the 
number of events concerned was not sufficient to affect 
the standard deviation noticeably. 

Double-Compton events in which both vertices are 
within the active volume produce two primary signals 
which are time coincident, but separated in position, 
and two secondary signals which are separated in both 
time (delay) and position. Even if they cannot be sepa- 
rated they are of no consequence as the combined ratio 
of S2/S1 will be relatively unaffected. However, if one of 
the vertices occurs in a position from which no secondary 
is possible, then the only way to identify them is through 
positional mismatch between SI and S2 and a less well 
reconstructed position from SI as this has two vertices. If 
the 'dead' vertex is very close to one of the PMT surfaces 
the SI signal can also appear to be too peaked within the 
array. Although there were already specific software cuts 
designed to deal with these events, some with certain 
topologies were not being fully identified by our analy- 
sis at that stage. For the 137 Cs data this problem was 
most apparent in the region log 10 (S2/S1) < —0.5 and 
E > 30 keVee but extended right down to the lowest 
energies. The 137 Cs calibration data were thus used to 
improve our algorithms for identifying MSSIs and the 
new routines were implemented after the science data had 
been opened. However, even with the improved selection 
cuts it was still not possible to use the 137 Cs data to 
predict accurately the expected number of single-scatter 
events leaking into the nuclear recoil region. The com- 
bination of lowering of the band mean (due to the rate 
dependent PMT sensitivity suppression at low temper- 
ature) towards the nuclear recoil band and remaining 
additional events in the lower wing caused a large over- 
prediction of event leakage into the WIMP search box (41 
events were predicted). The additional events remaining 
in the lower wing were probably due to the 137 Cs source 
not accurately mimicking that of the background sources 
due to its location. Hence, instead of using the 137 Cs 



6 



2-3 keVee 



7-8 keVee 



15-16 keVee 



5 



10' 7 



10 7 



x =0.073±0.011 
a=0.216±0.009 
a=1.854±0.261 
A=27.801±0.798 




5 



10 7 




%ldof=18.5/9 
H=-0.231±0.007 
a=0.131±0.006 
A=66.49±4.71 



10* 7 



10 7 



1 



io- 



■0.13, a=0.18 

x =-0.303±0.009 
io=0.24S±0.009 
a^l.S27±0.187 
A=69.829±1.648 



10 



-1 -0.5 0.5 

logJSZ/Sl) 




-0.5 0.5 

logjsmi) 



0.5 0.5 

log 10 (S2/Sl) 



FIG. 6: Statistical fitting of the electron and nuclear recoil populations using the WIMP-search (upper panel) and AmBe data- 
sets (lower panel). Three 1-keVee wide bins are shown: lowest, intermediate and highest energies accepted. The electron recoil 
population was fitted with a skew-Gaussian function using both the minimum \ 2 (thin blue line) and a maximum likelihood 
(ML) method with the Poisson distribution as estimator (thick red line). The latter fit is more appropriate to data with low 
statistics (including zeros) in the tails of the populations as it uses a Poisson distribution as an estimator. Note that the entire 
population can be fitted (all energy bins, across the entire log 10 (S2/S1) range). The ML best fit parameters are indicated, 
along with the mean and standard deviation of the skew-Gaussian. The lower panels show the log-normal fits to the AmBe 
recoil data, which is used to define the acceptance region [fi-2a,fi], between the vertical dashed lines. The number of electron 
recoils observed to be leaking into this region, n b 3 , is compared with the estimated number, n ca i, from the ML fits. The total 
number of events expected in the acceptance region is 11.6±3.0. 



data, the WIMP-search data themselves were used to es- 
timate the expected electron- recoil backgrounds, and this 
gave 11.6±3.0. 



III. ANALYSIS OF THE WIMP SEARCH DATA 

A. Data processing and selection 

The raw data were reduced using the purpose- 
developed code ZE3RA (ZEPLIN-III Reduction and 
Analysis). The DAQ hardware records the 62 wave- 
forms at 500 MS/s (2 ns samples) for 36 fis periods. 
ZE3RA finds candidate pulses in individual waveforms 



by searching for 3 rms excursions above the baseline. 
Subsequent waveform processing includes resolving adja- 
cent/overlapping pulses and grouping of statistically con- 
sistent structures (e.g. scintillation tails). A statistically- 
motivated timing/shape coincidence analysis was then 
used to correlate occurrences on different channels thus 
allowing further pulse interpretation (e.g. clustering, 
identification of random coincidences, etc.) The result- 
ing pulses were ordered by decreasing area in the high- 
sensitivity (HS) sum channel and the largest 10 were 
stored in databases for further analysis. By design, 
ZE3RA does not ascribe physical meaning to pulses, 
it rather parameterises them in terms of arrival time, 
width, area, amplitude, etc. An event browser allows vi- 
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FIG. 7: Double-logarithmic plot of Fig. 5 showing the nu- 
clear recoil population obeying the power-law trend indicated 
by the yellow line; the behaviour of the inelastic line from 
129 Xe is markedly different, as this is dominated by charge 
recombination of the 40 keV 7-ray rather than the small nu- 
clear recoil component of the deposited energy. Approximate 
thresholds for SI (3-fold software trigger) and for S2 (hard- 
ware trigger) are also indicated. From this it can be seen 
that the S2 hardware trigger corresponds to Sl=0.5keVee for 
nuclear recoils. 
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sual scanning of events, channels or individual pulses; a 
batch-mode interface allows scripted reduction of large 
data-sets. 

The data structures produced by ZE3RA were anal- 
ysed by a flexible code based on hbook [13]. ft pro- 
cessed the original parameters to assign physical mean- 
ing to pulses in events according to a well defined set 
of rules (e.g. primary scintillation signals are fast and 
must precede wider electroluminescence signals). Only 
events that can represent single scatters in the two-phase 
target ('golden' events with one SI and one S2) were re- 
tained. Primary (SI) pulses were found by applying an 
acceptance threshold of 1/3 p.e. to the ZE3RA pulses 
and also requiring a 3-fold coincidence amongst the 31 
PMTs. This software threshold was nominally equivalent 
to an energy threshold of 1.7keVee. Exceptions in the SI 
selection were allowed for PMT after-pulses. These are 
signal-induced artifacts generated within the PMTs. In 
general they have a characteristic time delay from the op- 
tical signal, but with a wide distribution and, moreover, 
it varies between PMTs. As a result it is not trivial to 
identify after-pulsing and avoid them instead being clas- 
sified as additional SI signals, which would result in the 
event being wrongly rejected. Secondary (S2) pulses were 
required to have at least an integrated area corresponding 
to the signal expected from about 5 electrons leaving the 
liquid surface. This suppresses optically-induced single- 
electron emission [22] as well as optical feedback effects 
from the cathode grid, which are not part of the direct 
measure of the ionisation signal generated at the inter- 
action site. Many additional parameters are derived for 
these, such as 3-D position information, hit-pattern de- 
scriptors, interaction energy and corrections (e.g. array 
flat-fielding, electron lifetime, liquid level, light collec- 
tion, etc). Subsequent analysis (science exploitation) is 
based on PAW [14] and ROOT [15]. 
Trapping MSSI events effectively was a significant chal- 
lenge, involving a combination of approaches: use of 
goodness of fit indicators in the position reconstruction 
algorithms, comparison of coordinates derived indepen- 
dently from SI and S2, and searching for abnormal light 
patterns across the array. 



B. The WIMP Search Box 

Discrimination between nuclear and electron recoils 
is illustrated in Figure 9 which combines electron re- 
coil data from 137 Cs and elastic nuclear recoil data from 
AmBe. The separation between the two populations is 
clear and this is used as the main way of defining the nu- 
clear recoil search box for potential WIMP events. The 
selection cuts used can be categorised as follows: 



FIG. 8: Comparison of the skew-Gaussian mean and standard 
deviations for the 137 Cs and WIMP-search data-sets calcu- 
lated from the ML-fit parameters (the horizontal error bars 
indicate the bin width). 



1. Golden event selection (including pulse finding, SI 
and S2 definition, and single scatter selection) 

2. Waveform quality cuts (mild cuts mainly aimed at 
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FIG. 9: Combined scatter plot of log 10 (S'2/Sl) as a function 
of energy from the two calibration data-sets, 137 Cs and AmBe. 
The upper population corresponds to low-energy Compton 
electrons and the narrower, lower one to nuclear recoils pro- 
duced by neutron elastic scattering. 



large baseline excursions compromising pulse pa- 
ramcterisation) 

3. Pulse quality cuts (mild cuts to avoid extreme out- 
liers in parameter distributions) 

4. Fiducial volume definition (drift time window and 
a radial limit from the S2 position reconstruction) 

5. Event quality cuts (strong cuts to deal with MSSI 
events mainly) 

The fiducial definitions (4) leave an active mass of 6.52 kg 
with a raw exposure of 453.6 kg. days. Low-energy events 
in the 10% data were well separated from the nuclear re- 
coil median line down to the lowest energies. The WIMP 
search box boundary was thus defined as 2<_E<16 keVee 
and (fi n — 2a)<log 10 (S2/Sl)<ii n , where /j, n is the energy- 
dependent mean of the nuclear recoils (acceptance of 
47.7%). This region was defined before unblinding and 
was kept for the subsequent analysis. The effective total 
exposure within this box, after taking account of all of 
the efficiencies, as detailed in Table I, is 127.8 kg-days. 



C. Backgrounds 

Electron and nuclear recoil background predictions for 
ZEPLIN-III are based on a full GEANT4 [16] simula- 
tion including measured radioactive content levels for all 
major components [3]. The largest contributor, by far, is 



FIG. 10: Electron recoil background measured during the 
fully-shielded science run. The differential spectrum is 
shown superimposed on the Monte Carlo prediction [3] us- 
ing GEANT4 [16] without rescaling. The latter includes a 
dominant 10.5 evts/kg/day/keVee ('dru') from the photomul- 
tipliers, 7-rays from the lead 'castle' (0.7 dru), /3-particles 
from 85 Kr (0.2 dru) and 7-rays from ceramic feedthroughs 
(0.1 dru). The disagreement at high energies is caused by 
single-scatter selection in the data (but not in the simulation) 
and by the limited DAQ dynamic range which was optimised 
for the WIMP-search run. 



the PMT array. Figure 10 shows the measured differ- 
ential background spectrum together with the simulated 
background. The high-energy region above 300 keVee is 
suppressed due to dynamic range limitation. 

The expected single-scatter neutron background in the 
data-set is 1.2 ± 0.6 in the WIMP search box with 90% 
coming from PMT generated events through (cv,n) inter- 
actions and spontaneous fission. The remaining 10% are 
mainly from contaminants in ceramic feedthroughs and 
external leakage through the shield of neutrons from the 
rock. 



D. WIMP Signal Search 

Figure 11 shows the final scatter plot from the com- 
plete science data-set. There are 7 events within the 
WIMP search box and the energy scale is shown in keVee. 
To assess the implications of these events the energy scale 
needs to be converted into keVnr, the energy dependent 
detector efficiency for nuclear recoils must be found and 
the relative likelihood of any of those 7 events being 
drawn from the expected WIMP distribution rather than 
the extended electron-recoil distribution must be calcu- 
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FIG. 11: Scatter plot of log 10 (52/51) as a function of energy 
for the entire 83-day data-set of first science run. There are 
7 events (large dots) in the WIMP-search region (thick red 
box), which extends from 2 < E < 16 keVee and fj, n — 2a n < 
log 10 (52/51) < fi n , where fj> n is the energy-dependent mean 
of the nuclear recoils (thin red line bordered by the blue curves 
at ±lcr„). These are all located near the upper boundary, 
between ~5-15 keVee. 
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FIG. 12: Horizontal distribution of events in the energy range 
2-16 keVee for the science run data-set produced in the same 
way as those in Figure 3. The reconstructed location of the 
7 events in the acceptance region is indicated. The measured 
distribution of the overall background is consistent with de- 
tailed Monte Carlo simulations of that expected from the in- 
strument activity which is dominated by the photomultipliers. 



lated. 

The level of discrimination apparent in Figure 11 is 
very high. As derived from the data themselves, the av- 
erage 7-ray rejection factor is 5xl0 3 between 2-16 keVee 
with an increase below 5 keVee. This is significantly 
better than had previously been demonstrated by the 
XENON10 experiment which achieved 99.9% at the very 
lowest energies [10] whilst our data exhibit better than 
99.99% in the 2-5keVee band. 

Figure 12 shows the spatial x-y distribution of all 
events in the 2-16 keVee energy range. Events within 
the WIMP search box are highlighted. 

To derive the significance of the events within the 
search box the experiment efficiency must be derived to- 
gether with the energy scale conversion between keVee 
and keVnr. These are established in the following sec- 
tions. First the efficiency for nuclear recoil detections is 
found by comparing AmBe data-sets with very different 
trigger thresholds in both hardware and software. The 
energy scale conversion is then done by comparing a sim- 
ulation of the expected nuclear recoil with that measured. 



1. Efficiency and threshold 

The overall detection efficiency will be a combination 
of hardware and software effects. As mentioned earlier 



in Section IB the hardware trigger threshold is derived 
from S2 in the low-energy part of the SI spectrum rele- 
vant to WIMP signals. At higher energies, well beyond 
the upper limit of the WIMP search box, there is a high- 
level inhibit to suppress the overall count rate, but this 
does not affect the efficiency at low energy. Dead-time 
effects are usually energy independent. Software effects 
include thresholding associated with pulse finding algo- 
rithms and selection cuts. These have been described 
in Section IC. The energies (expressed in SI keVee) 
at which these 'thresholds' will affect the detection effi- 
ciency are tabulated in Table I. In order to confirm these 
expectations a second AmBe data-set was analysed as a 
check on the energy dependence near the threshold. This 
data-set had been acquired with a lower hardware S2 
trigger threshold. In addition, the 3-fold SI coincidence 
requirement was changed to 2-fold in this particular anal- 
ysis and all quality cuts removed or significantly relaxed. 
The overall effect of these two changes is shown in Fig- 
ure 13 by comparing the black histogram labelled 'Am-Be 
low-threshold data' with the blue shaded histogram la- 
belled 'Am-Be calibration data'. The difference between 
these two histograms is only noticeable below SI ^4 keVee 
as expected. A study of the smallest S2 events triggering 
the system in each run has shown directly that the trigger 
level in the two runs was ~11 and ^4 ionisation electrons, 
respectively. These numbers were calibrated against the 
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measured single electron spectrum for ZEPLIN-III fol- 
lowing the method already used for ZEPLIN-II [22]. The 
experiment efficiency during the science run is taken as 
the ratio of the two AmBe data-sets, shown in Figure 14. 

The full red curve labelled 'Simulation 
(Eee/Enr=2.09)' in Figure 13 shows a Monte Carlo 
simulation of the expected differential spectrum which 
should have been seen by the experiment assuming a 
constant ratio between SI keVee and SI keVnr. This 
simulated curve has not been corrected for instrument 
efficiency but even so it is clear that there is a departure 
from the experimental data below SI ~ 20 keV. Given 
that this mismatch extends so far in energy above 
any reasonable thresholding effects, it is interpreted as 
evidence for a non-linear scale conversion. 



2. Energy conversion 

A comparison between the differential spectrum seen 
during the nuclear-recoil calibration, using AmBe, and 
a Monte Carlo simulation has been used to derive the 
energy scale conversion between keVee and keVnr. This 
relies on the integrity of the simulation using GEANT4, 
which is very well established in general for elastic scat- 
tering of neutrons, and which has been further exten- 
sively validated as part of this work. Systematic effects 
related to the simulation of the experimental calibration 
were assessed. These included, amongst others: varia- 
tions in neutron source spectrum and source location in- 
side the shield; the effect of intervening and surrounding 
materials; simulation event selection; energy resolution 
smearing; coincident AmBe 7-rays; treatment of inelastic 
scattering in xenon. The Monte Carlo result at low recoil 
energies was very resilient to sensible variation of these 
parameters. A different Monte Carlo code [17] confirmed 
these results independently. Naturally, incorrect angu- 
lar cross-sections for elastic scattering off xenon could 
be invoked to explain the low-energy result, since en- 
hancing forward scattering would soften the recoil spec- 
trum. However, dedicated simulations confirmed the cor- 
rect implementation of the ENDF/B-VI evaluated data 
libraries [18] which underpin the GEANT4 low-energy 
neutron transport models. Both angular and energy- 
differential cross-sections were found to be in agreement 
with ENDF/B-VI and similar databases. An implemen- 
tation in GEANT4 of the more recent ENDF/B-VII data 
for xenon by the XENON10 team [12], aimed at explor- 
ing the causes of a similar effect observed by that exper- 
iment, found only minor differences in the recoil spec- 
trum produced by a similar neutron source. We have 
independently confirmed this conclusion. The compari- 
son between simulation and experiment for ZEPLIN-III 
is shown in Figure 13. The energy scale associated with 
the simulated data has been converted from keVnr to 
keVee in Figure 13 by simply dividing by 2.09, to allow 
for the combination of the relative nuclear-recoil scintil- 
lation efficiency to that of a 122 keV 7-ray at zero electric 
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FIG. 13: Differential energy spectra for the AmBe elas- 
tic recoil population in SI electron-equivalent units ( 57 Co- 
calibrated SI). The main calibration data (shaded blue his- 
togram) and the lower threshold data-set described in the text 
(black) are compared with the Monte Carlo simulation using a 
constant conversion factor between nuclear recoil and electron 
equivalent energies (solid red curve) . The ratio between these 
two curves is interpreted as an energy-dependent efficiency 
factor and occurs in the low-energy region where threshold- 
ing effects are expected. The dashed red curve is the result 
of the nuclear recoil non-linearity analysis described in the 
text, which results in the energy conversion indicated by the 
markers at the top of the figure. 



field, L e ff, and a suppression factor, S, which allows for 
the field-dependent variation in the scintillation output. 
These are used in the following equation: 



SI S e 

Ly L e ffS n 



(1) 



where S e and S n are the suppression factors in the scin- 
tillation output for 122 keV 7-rays and nuclear recoils, 
respectively, at the experiment operating fields. Note 
that in this equation the ratio Sl/L y defines the keVee 
unit. Above E nr ~20 keV the available experimental data 
for L e ff suggests it is constant at ~0.19 [19-21]. How- 
ever a variation in L e f f at low energy has been invoked 
to explain XENON10 neutron calibration data [12] and 
hence is allowed to vary in this work. 

In general the conversion between an electron- 
equivalent energy scale, in keVee, and a nuclear re- 
coil energy scale, in keVnr, is not necessarily linear 
and any non-linearity can be expressed mathematically 
through energy dependency in L e ff and/or S e /S n . Above 
E nr ~20 keV the available experimental data for L e // sug- 
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TABLE I: Energy-independent efficiency factors and thresh- 
olds due to hardware and software actions. Efficiency fig- 
ures are constant over the WIMP recoil range. Numbers 
following the entries refer back to the list of software oper- 
ations itemised in Section IIIB. The total effective exposure 
is 127.8 kg-days. 



Effect 


Efficiency 


Method 


Deadtime 


91.7% 


Measured 


Hardware upper threshold 


100% 


On-off compare 


ZE3RA pulse finding (1) 


96.0% 


Visual inspection 






Hand calculation 


Event reconstruction (2,3) 


91.9% 


Visual inspection 


Selection cuts (5) 


73.0% 


On-off compare 


WIMP box acceptance 


47.7% 


Calculation 


p.- 1 J' ( 

Hj fleet 


1 nresnold 


Ivletnoa 


Hardware (S2) trigger 


Sl= 0.5 keVee 


Two data-sets 






Visual inspection 






Modeling 






Pulser tests 


Software S2 area 


SKI keVee 


Calculation 






Scatter plots 


Software SI 3- fold 


Sl= 1.7 keVee 


Calculation 






Two data-set 






analyses 



"All thresholds are quoted here in terms of the SI signal in keVee 
for nuclear recoils. The equivalent nuclear recoil energy, keVnr, 
depends on the conversion between keVee and keVnr. For the re- 
lationship shown in Section HID 2, 11 ionisation electrons corre- 
sponds to < 7keVnr 
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FIG. 15: The derived energy-dependent behaviour of L e ffS n - 
The thick curve shows the best fit to the data, but other 
curves producing very similar goodness-of-fit indicators are 
obtained within the envelope shown. The constraints become 
very weak outside the energy range shown. 

gests it is constant at ^0.19. At lower energies the situ- 
ation is much less clear [12]. For S n there are no data on 
the energy dependence but rather there is a single value 
based on a measurement at 56 keVnr using a neutron 
beam [20]. This gives S n — 0.90 at our field and it is 
commonly assumed to remain constant over the whole 
energy range of WIMP nuclear recoils. If L e // and/or 
S n are not constant below ~20 keVnr this will cause a 
non-linearity in the nuclear recoil energy scale. 

In the following it is assumed that such non-linearities 
are responsible for the mismatch seen in Figure 13. The 
approach used is similar to that applied to the XENON10 
data [12]. Using a maximum- likelihood technique we 
have derived a non-linearity function which best matches 
the AmBe simulation to our neutron calibration spec- 
trum above ^2 keVee. The outcome of this process is 
shown as the dashed red curve in Figure 13. Figure 15 
expresses the nonlinearity in terms of the combined effect 
of L e ff and S n , with the latter referenced to 0.90. In Fig- 
ures 13 and 14 the top horizontal axes show the energy 
scale in keVnr to be compared with keVee on the bottom 
scale. The WIMP search box boundaries then translate 
to 10.7 and 30.2 keVnr. One consequence of the required 
non-linearity is a marked reduction in efficiency for nu- 
clear recoil detection below 15 keVnr. 



FIG. 14: Energy-dependent part of the nuclear recoil de- 
tection efficiency as deduced empirically by comparing the 
two experimental AmBe spectra shown in Fig. 13. The 'low- 
threshold' run was taken with a lower hardware trigger thresh- 
old; in addition, software quality cuts were relaxed, along with 
the SI 3-fold requirement. A fit to the data is shown, with 
the WIMP acceptance box indicated by the thicker portion of 
the line. The S2 hardware trigger in the low-threshold AmBe 
run was half of that used in the science data-set and thus cor- 
responded to SI = 0.25keVee, and the SI pulse finding algo- 
rithms only required a 2-fold detection above 1/3 p.e. giving 
a nominal software threshold of SI = 1.1 keVee. Hence above 
SI = 2 keVee (the lower WIMP acceptance box boundary) the 
'low-threshold' data-set has near-unity efficiency. 



3. Limit analysis 

The event box contains a large empty region with a 
small number of events close to where a tail from the elec- 
tron recoil distribution is expected. However, although 
there is a good fit of a skew-Gaussian distribution to the 
electron-recoil band above the WIMP search box, there 
remains systematic uncertainty about an extrapolation of 
this being used as an accurate estimator of the number 
of expected background events in the box. The fact that 
the best fit expectation exceeds the measured number of 
events might result in an artificially lower upper limit, 
as pointed out in [23]. This compromises any straight- 
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forward use of maximum-likelihood techniques and even 
the commonly- used Feldman- Cousins (FC) analysis [23]. 
Hence, a simpler, more transparent and conservative ap- 
proach is adopted based on a minimum of three pieces of 
information about the data. 

The first is the reasonable assumption that any ex- 
pected electron-recoil background will fall in the top 
part of the WIMP search box. Based on this assump- 
tion the box is divided into two regions which have sig- 
nificantly different probabilities of having electron-recoil 
background within them. This is done in Figure 16 after 
transforming the WIMP search box so that the verti- 
cal axis has a linear scale in nuclear recoil acceptance 
percentiles as derived from the AmBe calibration data. 
In this representation any WIMP nuclear recoil signal 
should populate the box uniformly, whereas the density 
of the electron recoil background is expected to decrease 
monotonically down from the top. A horizontal dashed 
line is shown which divides the WIMP search box into 
two regions such that the top area contains all the events. 
In the following analysis the fractional area in the lower 
region is denoted by /. 

The second is the observation that no WIMP event is 
seen in the lower region (n; = 0). 

Finally, it is possible that there may be up to 7 WIMP 
events in the upper region (n u ^ 7). 

A classical 90% one-sided upper limit for the WIMP 
expectation value in the whole box, /j,, is the value under 
which 10% of repeated experiments would return zero 
events in the lower box and up to 7 in the upper box. 
This is expressed in terms of Poisson probabilities as 

P(ni =0,n u ^ 7\fi) = 

7 

P(m = 0|/ M ) * p ( n » = - = 0J ( 2 ) 

i=0 

Over the range of values of / between 0.75 and 0.84 the 
calculated result is /j, = 2.30/ /. / = 0.84 is the maximum 
area allowed which just excludes all of the events. 

It turns out that, for the value of fi resulting from this 
calculation, the second factor in equation (2) is very close 
to unity regardless of the area fraction, /. This reflects 
the fact that the upper limit is driven almost entirely by 
the presence of the empty region and the value 2.30 is 
then recognised as the classical 90% upper limit on zero. 
It is then reasonable to assume that the two-sided 90% 
confidence interval for this particular data-set will also 
be driven by the empty box. In this case the upper limit 
to this interval will be at \i = 2.44/ f, with 2.44 being the 
corresponding 2-sided FC upper limit on zero [23] . Figure 
16 shows a dividing line with / = 0.8, which is adopted 
as a conservative boundary placement beyond which no 
background is likely. The 90% confidence interval upper 
limit is then fi = 3.05. With this extreme value of fj, 
there is a 54% probability that there are indeed no WIMP 
events in the upper region, a 33% chance of there being 
1 WIMP event and a 13% chance of ^ 2 WIMP events. 
The fact that the most likely scenario is no WIMPs in 
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FIG. 16: The WIMP search box with the vertical axis 
remapped onto nuclear recoil percentiles. This is done us- 
ing the S2/S1 distribution from the AmBe calibration data. 
The positions of the 7 events falling within the box are shown 
as well as other events just outside the box. The horizontal 
dashed line separates the box into two regions with an area 
ratio of 1:4. 

the data-set even with fj, = 3.05 implies that ji = is 
included within the 90% two-sided interval as the null 
event hypothesis becomes more and more likely as /U is 
reduced. 

The upper limit of 3.05 events is used to derive the 
upper limit to the WIMP-nucleon spin-independent elas- 
tic scattering cross-section as a function of WIMP mass. 
The signal energy distribution is obtained from the 
theoretical WIMP recoil spectrum [24], derived using 
the standard spherical isothermal Galactic halo model 
(pdm=0-3 GeVcm -3 , v o ~220 km/s, v esc =600 km/s and 
VEarth=232 km/s), detector response efficiencies and en- 
ergy resolution. The form factor is taken from [25]. The 
expected distribution in S2/S1 is determined from the 
neutron calibration. 

The final result for the 90% confidence interval upper 
limit to the cross-section, shown in Figure 17, has a mini- 
mum of 8.1 x 10~ 8 pb for a WIMP mass of 60 GeV/c 2 . In 
the mass range beyond 100 GeVc -2 this result comple- 
ments the XENON10 result and further constrains the 
favoured SUSY parameter space [26] from xenon-based 
experiments. Spin-dependent limits are presented sepa- 
rately [29]. 



IV. CONCLUSIONS 

An analysis of 847 kg-days of data from the first sci- 
ence run of ZEPLIN-III has resulted in a signal lower 
limit consistent with zero, and an upper limit on the 
spin-independent WIMP-nucleon elastic scattering cross- 
section of 8.1 x 10 _8 pb, at 90% confidence level. In reach- 
ing this result it was necessary to confront an unexpected 
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FIG. 17: 90% confidence interval upper limit to the WIMP- 
nucleon elastic scattering cross-section as derived from the 
first science run of ZEPLIN-III for a spin-independent in- 
teraction. For comparison, the experimental results from 
XENON10 [10, 27] and CDMS-II [28] are also shown. Note 
that the XENON10 curve is a 1-sided limit, corresponding ap- 
proximately to an 85% confidence 2-sided limit [10]. CDMS-II 
and our result are both 90% 2-sided limits. The hatched ar- 
eas show 68% and 95% confidence regions for the neutralino- 
proton scattering cross-section with flat priors as calculated 
in Constrained MSSM [30]. 

mismatch between the nuclear recoil spectrum shown in 
the AmBe calibration data and the Monte Carlo simula- 
tion. A careful and thorough analysis of efficiency fac- 
tors and threshold effects (including the use of alternative 
data-sets with different thresholds, systematic changes to 
software cuts and thresholds, visual scanning and manual 
analysis of large samples of data and modelling and di- 
rect verification of the performance of the DAQ) did not 
resolve this mismatch. As a more credible alternative 
explanation the possibility of a non-linearity in the nu- 
clear recoil energy scale has been studied. Non-linearity 
as such is not unexpected and, indeed it would be sur- 
prising if it did not exist at low energy, and a similar ap- 
proach has been used by others for xenon [12]. Using this 
analysis it has been possible to reconcile the data with 
a non-linearity setting in at the same energy as in [12] 
but with a more significant effect at lower energies. In 
itself this may not be surprising given the very different 



operating conditions within ZEPLIN-III and XENON10: 
the most obvious being that the electric field in the liquid 
is 6 times stronger in the former. Indeed, there are other 
clear differences in the performances of the two instru- 
ments. However, it is clear that the physics underlying 
the low-energy performance is poorly understood. This 
is true of both the response to electron recoils [11] and 
to nuclear recoils [12]. As a point of reference, if the mis- 
match between the AmBe simulation and the data were 
interpreted solely as an instrument efficiency, the effect 
on the upper limit would not have been dramatic (<40% 
increase) as this approach has a better effective threshold 
for nuclear recoils but a poorer efficiency. 

The analysis presented is not blind as one of the analy- 
sis routines was changed after opening of the full data-set 
as was the limit setting procedure. In applying the limit 
analysis no use was made of any background estimates 
(neither electron-recoil or neutron scattering) and this 
was done deliberately to avoid underestimating the up- 
per limit. 
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